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I.  INTRODUCTION 

Autoregressive  spectral  estimation  has  become  an  important  method  of 
spectral  estimation  in  recent  years  (see  Akaike  (1969),  Parzen  (197*0, 

Ulrycb  and  Bishop  (1975),  and  Beamish  and  Priestley  (1981),  for  example) 
despite  1)  a  continuing  discussion  of  the  problems  of  order  determination 
and  2)  a  lack  of  easily  applied  procedures  for  determining  confidence 
intervals  or  bands  on  the  function  being  estimated. 

In  this  paper  we  derive  asymptotic  100  (1-a)^  simultaneous  confidence 
bands  for  an  autoregressive  spectral  density  assuming  one  has  data  from  a 
finite  order  autoregressive  process  with  known  order  or  conditional  upon 
having  correctly  determined  the  order  if  it  is  unknown. 

The  bands  are  derived  in  section  2  and  their  implementat i on  described 
on  simulated  data  in  section  3« 

2.  SIMULTANEOUS  CONFIDENCE  BANDS 

The  basic  property  of  autoregressive  processes  that  we  shall  use  is 
that  the  reciprocal  of  the  autoregressive  spectral  estimator  is  a  linear 
combination  of  a  finite  number  of  asymptotically  normal  random  variables 
whose  asymptotic  covariance  matrix  is  easily  consistently  estimated.  Thus 
one  can  use  a  Scheff4  (1953)  type  projection  argument  to  determine  asymptotic 
confidence  bands  on  the  reciprocal  of  the  autoregressive  spectral  density 
(and  thus  on  the  density  itself)  at  all  frequencies. 

Let  Y  be  an  autoregressive  process  of  order  p  with  coefficients 
a(l) , . . . ,a(p)  and  noise  variance  o2,  i .e. 
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P 

l  a(j)Y(t-J)  -  e (t)  ,  t-0,+1,+2,... 

J-o 


where  a(0)  ■  1,  e  is  a  white  noise  series  of  uncorrelated  zero  mean  random 

variables  having  common  variance  a2,  and  the  zeros  of  the  complex  valued 
P  j 

polynomial  g(z)  *  £  a(j)zJ  are  all  greater  than  one  in  modulus.  Then  we 

j-0 

write  Y  <v  AR(p,a,o2)  and  note  that  the  covariance  function  R(v)  ■  E(Y(t)Y(t+v)) , 

1  °°  _  j  v 

v=0,+J,...,  the  spectral  density  function  f(<u)  =  yjf  I  R(v)e  u  ,  we  t“ir  , 

Vc“« 

and  the  parameters  o  and  a2  are  related  by 


f(ui)  =  — 2  *  o)e[-ir,if] ,  (2.1) 

g(e  ) 


l  ct(j ) R(j-v)  =  6  a2  ,  v>0,  (2.2) 

J-0 


where  is  the  Kronecker  delta.  If  on  the  other  hand,  Y  is  a  moving  average 
process  of  order  q  and  parameters  f}(l ) , . . .  ,8(q) ,  and  o2,  i  .e 


Y (t)  -  l  R(kMt-k)  ,  t-o+i,..., 

k=0 


where  8(0)  =  1 ,  we  have,  with  h(z)  *  J.  S(k)  z  * 

k=0 


R(v) 


q-lv| 

o 2  I  B(k)8(k+|v|)  ,  I v | <_  q 

k-0 


0 


I  v  |  >q 


1 


1  R(v)e"  1  Va) 
v=-q 


[R(0) 


+  2  l 


R(v)  Cos  vto] 


f-  |  h(elu))  |  2  ,  toe  [-it  ,it]  . 
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Thus,  if  Y^AR(p,a ,o2)  we  have  f(u)>0  for  all  w  and 


h(to)  =  fTwT*  =  ^2^'"  lg(e'W)l2  “  xT (u) y 


where  x^(to)  =  (■=!>  ~  Cos  to . 1  Cos  pto) ,  and  yT  =  (y(0),y(l) . y(p)), 

-  /IT  IT  IT 
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y(v)  =  —t~  l  a(j)ot(j+v)  ,  v=0 . . 

°  j-0 


Thus  the  reciprocal  of  f  i s  a  linear  combination  of  the  parameters 
y(o)  »y(0  » •  •  •  »y(p)  • 

Let  Y,v/AR(p,o,o2)  and  let  Y(l) , . . . ,Y(n)  be  a  sample  realization  from  Y. 
Let  0T  =  (aT,o2)  *  (o(l ) , • • . ,a(p) ,o2)  be  estimators  of  the  parameters 
0^  -  (ot^,o2)  *  (o(l) • ,a(p) ,o2)  that  have  the  same  asymptotic  distribution 
as  the  maximum  likelihood  estimators  of  0,  i .e  (see  Parzen  (1961)) 


Jn  (0-0)  Np+)  (O,$(0)) 


(2.3) 


where 


k 


♦  (e) 


o 

2  o4 


and  r  is  the  (pxp)  Toeplitz  covariance  matrix  of  Y ( 1 ) , . . . ,Y(p) .  Estimators 
P 

having  this  asymptotic  distribution  include 

(1)  The  Yule-Walker  estimators  (YWE) ,  i .e  the  solution  to  the 
Yule-Walker  equations  (2.2)  with 


R(v)  *|  J  Y (tl Y(t+v)  ,  v=0, . . . ,p 

t=l 

replacing  R(v)  • 

(2)  The  Least  Squares  estimators  (LSE),  i .  e  the  values  a(l ) , . . . ,a(p) 
minimizing 


T 

$ (a)  =  l  {Y(t)-a(l)Y(t-l)-...-a(p)Y(t-p)}2, 
t~p+l 

and  a2  ■  S(a)/(n-p)  . 

3)  The  Burg  estimators  (BGE),  (see  Ulrych  and  Bishop  (1975)  for  a 
description  of  these  estimators). 

Then  the  autoregressive  spectral  estimator  f  of  f  consists  of  evaluating 
(2.1)  with  6  replacing  6.  If  the  order  p  is  unknown  a  priori  there  are 
several  procedures  (see  Hannan  and  Quinn  (1979) >  for  example)  for  finding  a 
value  p  from  the  data  which  as  n-*»  has  probability  one  of  determining  the 
correct  order  p. 

Then  the  confidence  bands  are  given  by  the  following  Theorem: 


Theorem  2. 1  Let  Y ( 1 ) , . . . , Y ( n)  be  a  sample  realization  from  a  Gaussian 
AR(p,a,o2)  time  series  and  let  a(1 ) , . • • »a(p) »  and  a2  be  estimators  satisfying 
(2.3).  Let 
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Y (v)  -  l  o(j)a(j+v)  ,  v=0,...,p 

a2  j-0 

and  h(u)  =  — - —  =  [y(0)  +2  V  y(v)  Cos  vw] .  Let  8^  *  (a^,  1/a2), 

f(«)  2"  v-1 

a  t  t  ^ 

8  *  (a  ,  1/a2),  and  B(6)  be  the  (p+1)  *  (p+1)  matrix  having  (£.,m)  element 

Btm(§)  ■  <  o2yU-l)  »  m=p+l 

"~T  6,  (Jl,m)+^7—  60(i,,m)  ,  m*l . . 

a  i  a  L 

v 

i=l , . . . ,p+l ,  where  <5^  (fc,m)  =  a(m+(£-l))  if  m<p-(fc-l)  and  zero  otherwise, 
while  62(2, m)  =  a(m-(£-l))  if  m>fc-l  and  zero  otherwise.  Then: 

a)  /n  (y-y)  +  Np+)  (0,$  (y))  where  $(y)  =  B (8) $ (8) B ($T) ,  and 

t (b)  *  a2r’’  0 

0T  2/a4 

b)  Asymptotically,  as  n-*”,  the  probability  is  at  least  1-a  that 
simultaneously  for  all  u  e  [0,ir], 

- — - - <.  f(to)  <_  - — ^ - 

h(u>)  +  S(w)  h(w)  -  S(u>) 


where  if  h(u>)  -  S(u>)  <  0  we  use  infinity  as  the  upper  limit,  and 
X2 

S2(u>)  -  — xT(u>)f(y)x  U), 

with  x2  .1  being  the  upper  a  critical  value  of  a  Chi  square  distribution 
a » p"*  * 

having  p+l  degrees  of  freedom. 

Proof  To  prove  (a)  we  note  that  from  (2.3)  and  the  definition  of  6  we  have 
/n  ( 0—  0 )  S  Np+j  (0,  $  (e) ) .  Then  defining  c  -  I/a2  and  c  =  I/a2,  we  have  for 
v*0 , • • • , p , 


/n  [y  (v)  -y  (v)  ]  =  4tv2  £  /rT  {ca(j)a(j+v)  -  ca(j)a(j+v)} 

j~0 

=  4ir2  £  /n  (ca(j)a(j+v)  -  ca(j)a(j+v)  +  ca(j)a(j+v) 
j*0 


-  ca(j)a(j+v)  +ca(j)a(j+v)  -  ca(j)a(j+v)} 


p-v 

=  4it2  l  /n  (ca(j+v)  [o(j)-o(j)l  +  a(j)a(j+v)  (c-c) 
j=0 


+  ca(j)[a(j+v)  -  a(j+v)  ]  } 


Z  4 ir2  l  (ca(j+v)  /n  [a(j)-a(j)l  +  a(j)a(j+v)  /n  (c-c) 
j-0 


+  ca(j )  r/n  [a( j+v)  -  a(j+v)]} 


-  4*2Y  ca(j+v)  ^  ta(J)-a(j)] 
j-0 


+  4ir2/n  (c-c)  £  a(j)a(j+v)  +  4n2  ^  co(j-v)  /n  [o(J)-o(j)]t 

J-0  J-v 


-  ir  i  i  ft'-’  •  • 
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where  we  write  X  %  Y  to  indicate  that  the  sequences  {X  },  {Y  }  of  random 
n  n  n  n 

variables  converge  in  distribution  to  the  same  random  variables,  and  we  have 
used  Slutsky's  Theorem  (Rao  (1973),  p.  122)  repeatedly. 

Thus  we  have  shown  that  /n  (y~y)  ^  B(g)  (6~B)  and  (a)  follows. 

From  (a)  then  we  have 

n(ri)T  (h}  5  xp+i 

and  since  the  elements  of  £(y)  are  continuous  functions  of  the  elements  of 
Y  we  have  (Rao  (1973),  p.  388) 

n(y-Y)T  $  '  (y)(y-y)  +  Xp+j 

Thus  the  probability  is  (asymptotically)  1-a  that  the  true  parameter  y  lies 
inside  the  ellipsoid  defined  as  the  set  of  vectors  l  satisfying 

(y-0T  h(y-j^)  <_  l  , 

where  M  «  n  t  (y)/x2  and  X  i  s  the  upper  a  critical  value  of  a 

T  l  a,p+l  a,p+l 

Xp+i  distribution.  But  (Scheffe  (1959),  p.  ^07,  et^.  sec. ) ,  y  is  in  this 
ellipsoid  if  and  only  if 

|XT  (y-y)|  <  [xT 

for  all  (p+1)  dimensional  vectors  x  and  thus  In  particular  only  if  (for 
vectors  x(u)) 
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|h(w)-h(uj)  |  <_  S(o>)  ,  [0,7t], 

giving  a  100(1-a)i  simultaneous  confidence  band  for  h(u>)  *  l/f(u>)  as 
h (oj)  -  S  (co)  h(co)  <_  h ( u) )  +  S ( u>)  . 

Then  the  reciprocal  of  this  gives  a  simultaneous  confidence  band  for  f, 
where  if  h(co)  -  S(w)  <  0  we  use  infinity  as  the  upper  limit  on  f(w)  which 
we  can  do  without  decreasing  the  probability  content  of  the  bands. 

3.  EXAMPLES 

To  illustrate  the  method  of  section  2  we  calculated  the  autoregressive 
spectral  estimator  and  confidence  bands  for  10  sample  realizations  of 
length  200  for  each  of  three  AR  processes  (see  Beamish  and  Priestley  (1981)); 
J_Le  AR(2,  -.4,  -.45,  1),  AR(5,  1.7,  2.4,  1.634,  .872,  .168,  l),  and  AR(4, 
-2.7607,  3*8106,  '2.6535,  .9238,  1).  These  processes  have  characteristic 
polynomial  zeros  whose  moduli  are  given  by  (1.11,  2.00),  (1.12,  1.12,  1.19, 
1.19,  3.33),  and  (1.01998,  1.01998,  1.019798,  1.019798)  respectively,  and  have 
variances  2.66,  30.59,  and  761.3.  Further,  the  ratio  of  the  maximum  to 
the  minimum  values  of  the  spectral  densities  for  the  three  models  are 
approximately  111.5,  41,800,  and  6*106.  Thus  the  performance  of  our  method 
on  these  models  should  be  representative  of  its  performance  on  a  wide  class 
of  AR  models.  In  Figures  A-C  we  display  on  a  log  scale  typical  point 
estimators  and  35%  confidence  bands  for  the  spectral  density  of  the  AR(2), 
AR(5) ,  and  AR(4)  models  respectively  using  the  Burg  estimators.  The  true 
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spectral  density  (the  curve  with  the  x's)  is  superimposed  on  the  graphs. 

Each  curve  consists  of  connecting  values  of  the  various  functions  at  the 
lol  equally  spaced  points  between  0  and  n  (labelled  on  the  graph  0  to  1/2 
cycles  per  sampling  interval). 

On  the  graphs,  the  largest  noninfinite  value  of  upper  limits  is  assigned 
to  l/(h(u))  -  S  (ai) )  whenever  h(a))  -  S(oj)  <  0.  Examples  of  this  include 
frequencies  near  .12  for  the  AR(4)  model  and  those  around  .32  for  the  AR(5) . 
We  have  found  this  to  be  a  useful  diagnostic  in  practice  as  peaks  such 
as  those  in  these  simulated  data  may  in  fact  be  indicating  the  existence 
of  deterministic  components  in  real  data  in  which  case  the  spectral  density 
does  not  exist. 

In  summary  then,  the  confidence  bands  given  in  section  2  appear  to 


work  well  on  a  wide  variety  of  autoregressive  models. 
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FIGURE  B.  Typical  Point  Estimator  and  95%  C 
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